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Polymer solutions subject to pressure driven flow and in nanoscale slit pores are systematically 
investigated using the dissipative particle dynamics approach. We investigated the effect of molecular 
weight, polymer concentration and flow rate on the profiles across the channel of the fluid and 
polymer velocities, polymers density, and the three components of the polymers radius of gyration. 
We found that the mean streaming fluid velocity decreases as the polymer molecular weight or/and 
polymer concentration is increased, and that the deviation of the velocity profile from the parabolic 
profile is accentuated with increase in polymer molecular weight or concentration. We also found 
that the distribution of polymers conformation is highly anisotropic and non-uniform across the 
channel. The polymer density profile is also found to be non-uniform, exhibiting a local minimum 
in the center-plane followed by two symmetric peaks. We found a migration of the polymer chains 
either from or towards the walls. For relatively long chains, as compared to the thickness of the slit, 
a migration towards the walls is observed. However, for relatively short chains, a migration away 
from the walls is observed. 

PACS numbers: 47.61.-k, 47.50.-d, 61.25.Hq 



I. INTRODUCTION 

The transport of polymer solutions in small confin- 
ing geometries remains a long stranding topic with many 
unanswered questions. The understanding of the trans- 
port of polymer solutions is important to chromatogra- 
phy, electrophoresis, adhesion, lubrication, polymer pro- 
cessing, oil recovery, etc. Recently, there has been an 
increasing interest in the areas of macromolecular trans- 
port in channels with widths comparable to size of a 
single molecule Q. For example, the understanding of 
the structure and dynamics of polymeric molecules in 
micro-channels is important to DNA sequencing in chan- 
nels with widths ranging from 10 to 50/xm @, DNA de- 
livery through micro-capillaries in gene therapy, and to 
lab-on-chip applications that involve polymers Is- 
sues of particular interest pertinent to polymer solu- 
tions in the presence of laminar flow is the mass dis- 
tribution of polymers across the channel, the polymers 
conformational distribution, and the effect of the poly- 
mer chains on the profile of the solution velocity field 

Qaiiisiimmmmmm. 

Past experiments were performed to investigate the 
conformational changes of polymers in the presence of 
elongational flow by Perkins et al. Ilal and in the pres- 
ence of shear flow by Smith et al. [TtJ through tracking 
fluorescently labeled DNA molecules. These experiments 
showed that the polymer chains become elongated along 
the flow direction. Several experiments indicated the ex- 
istence of a depletion layer in the polymers near the con- 
fining walls Q, y], H, 0, ■ The recent experiment by Horn 
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et al. 6] on a dilute solution of a high molecular weight 
polyisobutene in a very viscous fluid observed that the 
polymer solution exhibits flow with an apparent slip at 
the boundary. 

It is well known that at equilibrium conditions, poly- 
mer chains near walls, and in the non-adsorbed regime, 
are sterically depleted from the wall due their reduced 
conformational entropy [Tsj . In the presence of a 
pressure-driven (Poiseuille) flow, the shear rate across 
the channel is non-uniform. As a result, one expects a 
higher concentration of polymers in the mid-section of 
the channel since, there, the local shear is lowest, lead- 
ing thereby to an apparent slip at the walls Q. Due to 
the additional stretching of the polymer chains near the 
walls, as a result of the shear stresses originating from 
finite shear rate, thermodynamic arguments predict that 
the thickness of the depletion layer is thinner when com- 
pared to that at equilibrium conditions. Consequently, 
the steric effect of the walls on the polymer chains is 
less severe than at equilibrium conditions. However, ex- 
periments on dilute polymer solutions demonstrated 
that even in the case of uniform shear flow, the polymer 
chains migrate away from the wall leading to a depletion 
layer that is wider than at equilibrium. Therefore, the 
thermodynamic arguments alone are unable to fully ac- 
count for the migration of the polymer chains observed in 
experiments. It is important to note that since these sys- 
tems are not at equilibrium, thermodynamic arguments 
are not applicable. 

The understanding of the migration of the polymer 
chains away from the wall must therefore take into ac- 
count hydrodynamic interactions in the solution. How- 
ever, such understanding is complicated by the non- 
Newtonian properties of polymer solutions. Theories 
have therefore been based on simplified systems rather 
than on polymer solutions. Previous theoretical and com- 
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putational studies that focused on the density profiles of 
polymers in uniform shear flow (i.e., Couette flow) or 
pressure driven flow f»',e,, Poiseuille flow) have been con- 
flicting n in m Il3l . Previous Brownian dynamics 
simulations showed a migration of the polymer chains 
towards the wall jjJEJ. In these studies, however, wall- 
monomer hydrodynamic interactions are not accounted 
for. A recent kinetic theory by Ma and Graham 12] of 
an infinitesimally dilute solution of dumbbells predicts a 
migration of the dumbbell particles away from the wall 
towards the mid-section of the channel. A key reason 
for the migration of the dumbbells away from the wall is 
attributed to the hydrodynamic interaction between the 
wall and the dumbbell, which leads to a deterministic 
migration from the wall, while the effect of the Brown- 
ian motion is to homogenize the solution and therefore 
to counter the effect of the wall-monomer hydrodynamic 
interaction. This cross-stream migration was recently 
investigated numerically by Jendrejack et at using a 
Brownian dynamics simulation that sclf-consistently ac- 
counts for both monomer- monomer hydrodynamic inter- 
action and wall-monomer hydrodynamic interaction . 
In this study, they observed that in the presence of pres- 
sure driven flow, the polymer density profile exhibits a 
local minimum (a dip) in the middle of the pore, which 
decreases as the flow rate is increased. They also ob- 
served that the depiction region from the walls is am- 
plified as the flow rate is increased. The depletion from 
the wall becomes substantial for large flow rates. In the 
more recent simulations by Usta et al. [TH . Ilif , which are 
based on a bead model for the polymer chains coupled to 
a lattice Boltzmann description for the solvent, the poly- 
mer density profile was also found to exhibit a minimum 
in the center-plane of the slit followed by two symmetric 
maxima. Near the wall, a migration away from the wall 
is observed when the ratio between the thickness of the 
slit and the equilibrium radius of gyration of the poly- 
mer chains is large, i.e., in the weak confinement regime. 
However, when confinement is increased, a migration to- 
wards the wall is observed. This result contrasts that of 
Jendrejack et al. |ll|. where even in the strong confine- 
ment regime, a migration away from the wall is observed. 
In the presence of uniform shear flow, Usta et al. [T^J ob- 
served that the density profile does not exhibit a peak in 
the mid-section of the slit, instead of a dip. The pres- 
ence of a dip in the mid-section of the slit must therefore 
be associated with gradients in the shear rate, present in 
pressure driven flow. 

Very recently, Khare et al. used molecular dynamics 
to study the density profile of dilute polymer solutions 
in both uniform shear flow and pressure driven flow |15| . 
In the case of pressure driven flow, they also observed a 
local minimum in the polymer density profile. However, 
the depletion layer is only weakly modified by flow. As in 
Usta et al.'s work [j^j. Kare et al. also observed that the 
polymer density profile peaks in the mid section of the 
slit in the case of uniform shear flow. A recent simulation 
by Fan et al. Q based on dissipative particle dynamics 



of dilute polymer solutions in pressure driven flow in a 
slit channel focused on the velocity profile of the poly- 
mer solutions. In their study however, they studied the 
migration of dimers instead of chains, and they observed 
very weak migration |8( ■ Therefore, the issue of migration 
clearly remains settled, warranting further investigation. 

A large number of studies during recent years have 
proved that dissipative particle dynamics (DPD) is a suc- 
cessful explicit computational approach in the investi- 
gation of both equilibrium and non-equilibrium proper- 
ties of complex fluids. A recent study by Ripoll et al. 
showed, in particular, that DPD correctly describe large 
scale and mesoscopic hydrodynamics in fluids [2(| • Us- 
ing extensive DPD simulations, we recently investigated 
equilibrium dynamics of polymer chains in the dilute and 
concentrated regimes, and both in bulk and in confining 
geometries. We found that the chains dynamics is well 
described by the Zimm model in the dilute regime, and 
as expected a crossover to Rouse dynamics is observed 
as the polymer concentration is increased [2ll . These 
studies validate the usefulness of DPD in correctly de- 
scribing the dynamics of polymer solutions. We there- 
fore used the DPD approach in studying the dynamics of 
polymer solutions in pressure driven flow. 

In the present article, we present results based on 
DPD simulations of polymer solutions in nanoscale slits 
and in the presence of pressure driven flow. Our study 
is performed on systems with various polymer volume 
fractions, chain length, and flow rate. We will particu- 
larly discuss the velocity profile of the solution across the 
channel, the conformational distribution of the polymer 
chains and their mass distribution. In Sec. 3, we present 
results of Poiseuille flow of a simple fluid using DPD. In 
Sec. 4, the model and computational approach is pre- 
sented. In Sec. 3, we present and discuss the numerical 
results. Finally, we summarize and conclude in Sec. 5. 



II. MODEL AND METHOD 

In this work, we investigate the dynamics of dilute 
and semi-dilute polymer solutions under pressure driven 
flow using the dissipative particle dynamics (DPD) ap- 
proach. DPD was first introduced by Hoogerbrugge and 
Koelman more than a decade ago |25j, and was cast in 
its present form about five years later [24l |25| . DPD 
is reminiscent of molecular dynamics, but is more ap- 
propriate for the investigation of generic properties of 
macromolecular systems. The use of soft repulsive inter- 
actions in DPD, allows for larger integration time incre- 
ments than in a typical molecular dynamics simulation 
using Lennard- Jones interactions. Thus time and length 
scales much larger than those in atomistic molecular dy- 
namics simulations can be probed by the DPD approach. 
Furthermore, DPD uses pairwise random and dissipative 
forces between neighboring particles, which are interre- 
lated through the fluctuation-dissipation theorem. The 
pairwise nature of these forces ensures the local conser- 



3 



a random force, . These three forces are respectively 
given by, 
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FIG. 1: (a) Velocity field along the a;-axis vs. z for the case 
of a driving force, f x — 0.02e/r c . Curves from bottom to top 
correspond to wall-solvent interaction parameter a ws = 3e/r c , 
10e/r c , 20e/r c , and 30e/r c , respectively, (b) The velocity 
profile vs. z 2 . Not that v x (z) is quadratic in z regardless 
of the wall-solvent interaction. The vertical dot-dashed lines 
indicate the positions of the walls. Notice that curves for 
a ws = 3e/r c and 10e/r c overlap with each other. 
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FIG. 2: (a) Velocity field along the a;-axis vs. z for a ws = 
3e/r c . Curves from bottom to top correspond to a driving 
force f x = 0.01, 0.02, 0.03, and 0.04e/r c . (b) The velocity 
profile normalized by the driving force vs. z. The excel- 
lent collapse in (b) implies that the velocity profile follows 
Eq. (11311 . (c) The velocity profile vs. z 2 . 
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where = r 3 - r h n i3 = r-y/ry, and = Vj - 
is a symmetric random variable satisfying 

<%(*)) = 0, 
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with i ^ j and k 7^ I. In Eq.©, At is the iteration time 
step. The weight factor is chosen as 



, > _ . 1 — r/r c for r < r c , 
^ ' ■ for r > r c 
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where r c is the interactions cutoff radius. The choice of 
uj in Eq. © ensures that the interactions are all soft and 
repulsive. The values of the parameter a%j used in the 
simulations are specifically chosen as 
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where w, s and p designate a wall, solvent and a monomer 
particle, respectively. We should not that a ww is irrele- 
vant since the wall particles are kept frozen in the simu- 
lation. The importance of the wall-solvent interaction to 
the velocity field boundary condition will be discussed in 
Section III. The integrity of a polymer chain is ensured 
via an additional harmonic interaction between consecu- 
tive monomers, 
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where we set, for the spring constant and the preferred 
bond length, C = lOOe and b = 0.45r c , respectively. 
The equations of motion of particle i are given by 
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vation of momentum, a necessary condition for correct 
long-range hydrodynamics [2(], . The DPD approach 
was recently used for the study of flow of colloidal sus- 
pensions [23 and flow of dilute polymer solutions [3, |28| . 

In the DPD approach, two particles, i and j, interact 
with each other via three pairwise forces corresponding to 



a conservative force, F^f \ a dissipative force, F\Y' , and 
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where m is the mass of a single DPD particle. Here, 
for simplicity, masses of all types of dpd particles are 
supposed to be equal. Assuming that the system is in 
a heat bath at a temperature T, the parameters 7^- and 
<Jij in Eqs. @ and J3J are related to each other by the 
fluctuation-dissipation theorem, 
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Consider a fluid in a pore of cross section A under a 
pressure gradient dP/dx. As a result, a portion of fluid of 
length dx along the x-axis experiences a net force along 
the direction of diminishing P given by, dF = AdPx. 
Consequently, each particle i in the volume Adx experi- 
ences an additional force 



f; = f x x = - - 
P 



dP 



dx 
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to be added to the forces in Eqs. (1-3). The values of the 
driving force considered in the present study range from 
to 0.04e/r c . 

Most of our simulations are performed on boxes of di- 
mensions L x x L y x L z with L x — L y — 19.6r c and 
L z = 20r c . We made sure that finite size effects are 
absent by performing simulations on boxes with larger 
values of L x and L y . As will be discussed in the next 
section, a frozen wall is used. The wall is parallel to the 
xy-plane, and the flow is along the x-axis. 

We used periodic boundary conditions in the three di- 
rections. Therefore only one wall is used. We consider 
a fluid with a number density p = 4r~ 3 , fc^T = e, and 
a = 3.0(e 3 m/r^) 1 / 4 . We used a wall with thickness 2.8r c . 
The wall is constructed from dpd particles that are ar- 
ranged in a face centered cubic lattice with number den- 
sity 61.35r~ 3 , and therefore much higher than that of 
the solvent. This is done in order to prevent the solvent 
particles from penetrating the wall. The iteration time 
5t = O.OIt, where the time scale r = (mr 2 /^) 1 / 2 . All our 
simulations were at least run over a time period of 8000t, 
corresponding to 8 x 10 5 time steps. The first lOOOr is 
used for relaxation of the system towards its steady state, 
and the remaining time is used for collection of the data. 



III. POISEUILLE FLOW OF A SIMPLE FLUID 
AND PARAMETER SELECTION FOR 
WALL-SOLVENT INTERACTION 

Before presenting results on the flow of polymer solu- 
tions, it is useful to investigate the flow of the solvent 
alone (a Newtonian fluid) in a slit, using the DPD ap- 
proach, and compare with the theoretical expectations 
for Poiscuillc flow. We recall that in a laminar flow along 
the x-axis, the xz-component of the stress tensor, H xz 
is related to the shear rate, (dv x /dz), through the rela- 
tion H xz — —rio(dv x /dz), where 770 is the fluid viscosity. 
The solution of the steady-state Navier-Stokes equation 
of the incompressible Newtonian fluid VP = 770 V 2 v, with 
a no-slip boundary condition at the walls, which are per- 
pendicular to the z-axis, v (z — ±H/2) — 0, and under 
a uniform pressure drop, (dP/dx) along the x-axis, cor- 
responds to the usual Hagen-Poiseuille parabolic profile 
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FIG. 3: Velocity field along the a;-axis vs. z for a WB = 3e/r c . 
Curves from bottom to top correspond to a driving force f x = 
0.01, 0.02, 0.03, and 0.04e/r c . Inset shows the same velocity 
profile normalized by the driving force. The excellent collapse 
of all data is an implication that the velocity profile follows 
Eq. fT3l . 



Within the DPD and molecular dynamics approaches, 
several techniques have been proposed to provide a no- 
slip wall- fluid boundary conditions. Revenga et ai, for 
example, used an effective field, instead of a wall com- 
posed of particles, in conjunction with a reflection force 
to reflect particles that cross the wall 30] . Willemsen et 
al. . [3l| proposed a method that involves a layer of par- 
ticles along the wall outside the simulation box. The co- 
ordinates and velocities of these particles arc determined 
by mirroring those of the fluid particles within the cut- 
off distance of the interaction potential. The momenta 
of the additional particles are such that the average ve- 
locity of a fluid particle, within the cutoff distance, and 
its auxilary particle within the wall satisfies the desired 
boundary condition. In a recent study, we found that 
the use of a thermalized wall can also help in obtaining a 
no-slip boundary condition |32j| . In the current study, we 
found that a no-slip boundary condition can practically 
be obtained by lowering the wall-fluid interaction. 

In Fig.QJ the velocity profile of a pure solvent is shown 
for the case of a slit width L z = 20r c , a driving force 
f x = 0.02e/r c and values of the wall-solvent interaction 
corresponding to a ws — 3e/r c , 10e/r c , 20e/r c , and 30e/r c . 
This figure shows that the velocity profile exhibits a dis- 
continuity when the amplitude is large (a ws > 10e/r c ), 
implying that the fluid lubricates the wall (i.e. the wall 
does not act as a no-slip boundary) at large wall-solvent 
interactions. However for a ws < 10e/r c , the wall acts ef- 
fectively as a no-slip boundary condition. In the remain- 
ing if this article, all simulations were performed with a 
wall-solvent interaction parameter a ws — 3e/r c . 

In Fig. Ufa), the velocity profile of the simple fluid is 
shown for different values of the driving force, and a wall 
solvent interaction, a ws = 3e/r c . As shown in Fig.^b), 
the normalization of the z-component of the velocity field 
by the driving force yields a universal velocity profile, in 
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FIG. 4: (a) Velocity field along the cc-axis vs. z for a polymer 
solution with N = 100 and at a driving force f x — 0.02e/r c . 
Curves from top to bottom correspond to polymer volume 
fraction ip p — 0, 0.06, 0.12 and 0.24, respectively. a ws — 
3e/r c . (b) Same data in (a) plotted vs. z 2 . The straight lines 
are linear fits to the data for small values of z. 



excellent agreement with the theoretical prediction for 
Poiseuille flow, Eq. (|T3|l . We also tested the scaling of 
the velocity profile with the width of the slit pore. As 
shown in Fig. |3J we found again that the scaling of the 
velocity profile with the width of of the slit is in excellent 
agreement with Eq. I|13|) . 



IV. POISEUILLE FLOW OF POLYMER 
SOLUTIONS 
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FIG. 5: Fluid velocity profile across the slit for the case of (a) 
•p = 0.12 and N = 50 and (b) ip = 0.24 and N = 100. In (a) 
and (b) Curves from bottom to top correspond to f x = 0.005, 
0.01, 0.02, and 0.04e/r c . The velocity profiles scaled by the 
driving force are shown for (c) p> — 0.12 and N = 50, and (d) 
ip — 0.24 and N — 100. The velocity profiles scaled by f^ n 
are shown in (e) <p = 0.12 and N = 50, and (f) <p — 0.24 and 
TV = 100. In (e), n = 0.92 and in (f), n = 0.84. The vertical 
dot-dashed lines indicate the positions of the walls. 



Velocity Profile of Polymer Solutions in 
Poiseuille Flow 



For a polymer solution, which is a non-Newtonian 
fluid, the shear stress relation to the shear rate, (dv x /dz), 
is approximated by a power law equation proposed by 
Ostwald [33| . 



solution, 



n*. = -k 



dv x 
dz 



(14) 



where K and the exponent n, are usually referred to as 
the flow consistency index and the fluid power-law index, 
respectively. For a polymer solution, which is a shear- 
thinning fluid, the power-law index n < 1. Eq. (1) trans- 
lates into a solution viscosity ry ~ (dv x / dz)"" 1 . Note 
that Eq. (1) implies that the viscosity diverges for very 
low shear rates. Therefore, the power law model is not 
expected to hold at very low shear rates. 

The balance between the driving force due to pressure 
gradient and the forces due to shear stresses lead to the 
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which reduces to Eq. (|13fl when n=l. 

We performed a large number of systematic simula- 
tions of flow of polymer solutions for different values 
of chain length corresponding to N = 5, 20, 50, and 
100, polymer volume fraction corresponding to tp = 0.06, 
0.12, and 0.24, and driving force ranging between and 
0.04e/r c . In Fig.0J the velocity profile of a polymer solu- 
tion with chain length N = 100 and under a driving force 
f x = 0.02e/V c is shown for volume fractions, ip = 0, 0.06, 
0.12 and 0.24. This figure shows that except from the re- 
gion close to the wall, the velocity profile is not quadratic 
in z, and therefore does not follow Hagen-Poiscuillc's law 
for Newtonian fluids, Eq. i|13|) . This is due to the non- 
linear rheological properties of the polymer solution. In 
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vicinity of the wall, however, the velocity profile is prac- 
tically independent of tp. This is presumably due to the 
fact that the polymer chains are sterically depleted from 
the wall regions. In Fig.^Jb), an approximate linear fit 
of the velocity profile vs. z 2 around the mid-section of 
the slit leads to an effective wall boundary that extends 
beyond the physical boundary. This leads to the an ap- 
parent slip of the polymer solution at the wall. 

In order to investigate even further the departure of 
the velocity profile from that Hagen-Poiseuille quadratic 
profile of Newtonian fluids, the velocity profiles of two 
polymer solutions, corresponding to (<p,N) = (0.12,50) 
and (0.24, 100) are shown in Fig. Fig. EJc) and (d) 
show that, in contrast to the case of a simple fluid (see 
Fig-Etb)), the velocity profile of a polymer solution does 
not scale linearly with the driving force, and that the 
departure from this scaling worsens as the chain length 
and/or the volume fraction of the polymers is increased. 
It is interesting to note that the scaling holds reasonably 
well close to the walls. This is due to the very low poly- 
mer density near the walls. A fit of the velocity profile 
with the generalized Hagen-Poiseuille equation for Non- 
Newtonian fluids, Eq. I|15|) . yields power law fluid index 
n = 0.92 and n = 0.84 for the case of (ip, N) = (0.12, 50) 
and (0.24,100), respectively. A modified scaling of the 
velocity profile with f^ n is depicted in Figs. |S[e) and 
5(f) respectively. Fig. [SJe) shows that the scaling holds 
reasonable well in the case of relatively low volume frac- 
tion of polymer. Fig. (Sff) shows that the scaling does 
not hold as well in the case of (ip, N) = (0.24, 100), where 
both the polymer volume fraction and molecular weight 
are higher. This is presumably due to the fact that the 
exponent n itself depends on the driving force, f x . 

In Fig. El the velocity profile for a solution with 
(<p,N) — (0.12,50) under a driving force f x = 0.02e/r c 
is shown for two values of the slit thickness correspond- 
ing to H = 7.8r c and 17.6r c . Fig. HJb) shows that the 
scaling of the velocity profile with H 2 in the case of a 
polymer solution does not hold in the center-plane re- 
gion, but holds reasonably well in vicinity of the walls. 
This is again due to the depletion of polymers from the 
walls region. Fig. |SJc) shows that a better scaling in the 
mid-section region is achieved when Eq. (|15fl is used with 
n = 0.84. 

We now turn our attention to the analysis of the poly- 
mers velocity profile. Since the polymers are sterically 
depleted from the regions near the wall, it is expected 
that in average the polymer chains move faster than the 
solvent and that the longer polymer chains move faster 
than the shorter polymer chains |3J| . In Fig. 7, the ratio 
between the average polymer velocity and the average 
fluid velocity, (v%)/{v x ) is shown. This figure shows that 
the polymer moves faster than the solution by about 10%. 
This effect is however reduced as the polymer volume 
fraction is increased. This is due to the fact that across 
the slit, the polymers become more evenly distributed 
as their volume fraction is increased. Fig. 7 shows as 
well that as the chain length is increased, (v^,)/(v x ) in- 
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FIG. 6: (a) Fluid velocity profile across the slit for the case 
of if = 0.12, N = 50 and f x = 0.02e/r c . Solid circles and 
open circles correspond to a slit thickness H = 7.8r c and 
H — 17.6r c , respectively, (b) The fluid velocity scaled by H 2 
vs. z/H. (c) The fluid velocity scaled by H^ +1 vs. z/H, 
with n = 0.84. 
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FIG. 7: The ratio between average polymer velocity and solu- 
tion velocity along the flow direction, (v x )/{v x ) as a function 
for chain length for f x — 0.02e/r c and H = 7.8r c . Curves 
from top to bottom correspond to <f> = 0.06, 0.12 and 0.24, 
respectively. 



creases, mainly as a result of exclusion of the polymers 
from the walls. However, for N = 100, Fig. 7 shows that 
(v%.)/ (v x ) is slightly smaller than that for shorter chains. 
Presumably, this is due to the fact that long chains tend 
to migrate away from the center-plane, as will be discuss 
later. 

In Fig. [SI the velocity profile of the polymer chains 
is shown for the case of N = 100 at a driving force 
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FIG. 8: (a) Polymer velocity field along the x-axis vs. z for 
a polymer solution with N = 100 and at a driving force f x = 
0.02e/r c . Curves from top to bottom correspond to polymer 
volume fraction <p p = 0.06, 0.12 and 0.24, respectively. a wa — 
3e/r c . (b) Same data in (a) plotted vs. z 2 . 

fx = 0.02e/r c for volume fractions p = 0.06, 0.12 and 
0.24. This figure shows that the velocity profile of the 
polymer chains differs from that of the embedding sol- 
vent despite the fact that the polymer chains are trans- 
ported by the solvent. Fig. EJa) and (b) depict the ve- 
locity profiles of the fluid and the polymer for the case of 
(<p,N) = (0.06,5) and (0.06,100), respectively. There, it 
is clear that the two profiles do not coincide, particularly 
when the chain length or volume fraction is increased. 
Fig. [5] shows that the polymer velocity exhibits a sharp 
discontinuity near the walls of the slit, implying that as 
far as the polymer is concerned, the walls act as a slip 
boundary. Fig. [5] shows in the middle of the slit, the 
polymer chains move slightly slower than the embedding 
solvent, and that the difference between the solvent veloc- 
ity and the polymer velocity is increased as the polymer 
chain length or volume faction is increased. Fig. [§] also 
shows that in the middle of the slit the polymer motion 
is retarded as compared to that of the solvent. However, 
the motion the polymers is faster than the solvent near 
the walls. The slower motion of the polymer chains near 
the walls is due to the fact that due to the monomers 
connectivity, monomers belonging to a single chain have 
to move with same average velocity. 



B. Polymers Conformational Distribution 

We now focus on the effect of flow on the conforma- 
tional distribution of the polymer chains in solution. In 
Fig. ^3 the distributions of the three components of the 
radius of gyration, Rg, R g v ^ and R g z \ of a polymer so- 
lution, with N = 50 and ip — 0.12 are shown for driving 
forces corresponding to f x = 0, 0.01e/r c , 0.02e/r c , and 




FIG. 9: Solvent velocity profile (dashed line) is shown to- 
gether with the polymer velocity profile (solid line) for the 
case of tp = 0.06 and f x = 0.04e/r c . (a) and (b) correspond 
to TV = 5 and TV = 100, respectively. The vertical dot-dashed 
lines indicate the positions of the walls. 



0.04e/r c . The z-axis represents the position of the center 
of mass of the chains. This figure indicates that, near the 
walls, the z-component of the radius of gyration is lower 
than than R g x ^ and R g v ^ near the walls. This is due to the 
steric interactions with the walls. Fig. 1101 shows that at 
equilibrium (f x — 0) and away from the walls, the three 
components of the radius of gyration are equal, implying 
that the chains away from the walls do not experience 
confinement. This implies that under non-equilibrium 
conditions, the redistribution of the radius of gyration 
is due to the flow. As the driving force is increased, 
the chains become stretched along the flow direction, i.e. 
along the ir-axis, regardless of the location of their cen- 
ter of mass across the channel. This figure also shows 
that the chains stretching along the flow direction is low- 
est in the mid-section of the slit pore. This is due to 
the fact that the shear rate, and therefore the resulting 
stresses, are lowest in this region. Away from the mid- 
section of the slit, the chains stretching increases but 
in a non-monotonous fashion. In particular, we observe 

(x) 

that the distribution of R g exhibits two peaks some- 
what close, but not very close, to the wall. The positions 
of these two peaks along the z-axis approach the wall as 
the driving force is increased. Very close to the walls, 
and within the depletion layer, the stretching is again 
increased. A similar complex behavior of the distribu- 
tion of chains radius of gyration was recently reported 
by Jendrejack et al. using their self-consistent Langevin 
dynamics approach We also note that the position 
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of the two weak off-center maxima in the distribution of 
the x-component correlate with the onset of the shoul- 
ders in the velocity profile of the polymer chains, shown 
in Fig.Efb). 

Fig. EI] shows that the y-component of the radius of 
gyration, R g , decreases from its maximum value at 
the center-plane up to some point, and then increases 
at distances closer to the wall. Since the shear stress 
11^ = 0, the contraction of the chains along the y-axis 
and around the center-plane is a partial compensation 
to the chains stretching along the x-axis. The increase 
of the y-component of the radius of gyration, closer to 
the walls, is a compensation to the decrease of the z- 
component of the radius of gyration close to the walls. 
As for the z-component, R g , its fast decrease near the 
walls is due to steric interactions with the walls, around 
the mid-plane, the decrease is slower, and is due to gradi- 
ent in shear rate along the x-axis. As the driving force is 
increased, shear stresses are amplified leading to further 
contraction of the chains along the z-axis. 

In Fig. ^2 t ne profiles of the components of the radius 
of gyration across the channel are shown for a solution 
with (p p = 0.12 with N = 5, 20, 50, and 100 at a driving 
force f x = 0.04e/r c . This figure shows that the effects 
observed in Fig. ^3 are enhanced as the chain length is 
increased. 



C. Polymers Density Profile 

In Fig. ll2f ah the monomers distribution in the case of 
a dilute solution with (p = 0.06 and N — 100 is shown for 
three values of the driving force corresponding to f x = 0, 
0.02 and 0.04e/r c . This figure shows that at equilibrium 
conditions (f x = 0) and away from the walls, the distri- 
bution is, as expected, essentially uniform. However, in 
the presence of flow, the distribution exhibits a local min- 
imum in the center-plane. This figure also shows that the 
thickness of the depletion layer slightly decreases as the 
driving force is increased. This figure therefore implies 
that for this specific tp and N, the polymer chains mi- 
grate away from the center-plane towards the walls. We 
note that in this case the ratio between the slit thickness 
and the chain radius of gyration is 2H/R g ra 3.6. For 
this ratio, Usta et al. also found a migration towards 
the wall [lj|- I n the case of weak confinement, <f> = 0.06 
and iV = 20 (for which 2H/R g ra 10.8), we also found 
a peak in the middle of the slit, although the peak is 
weaker than for N — 100. Near the walls, we found a 
slight increase in the depletion layer, again in qualitative 
agreement with Usta et aVs work [lj|- The presence of a 
dip, two off-center peaks and the slight deviation in the 
depletion layer near the walls was also observed in the 
recent molecular dynamics simulations by Khare et al. 

m 

In Fig. ll2f b). the monomers distributions for ip = 0.06, 
and chain length N = 5, 20, 50 and 100 are shown for a 
driving force f x = 0.02e/r c . This figure shows, that the 
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FIG. 10: Distribution of the three components of chains ra- 
dius of gyration along the slit for the case of iV = 50 and 
tp = 0.12. In (x), curves from bottom to top correspond to 
f x = 0, O.Ole/rv, 0.02e/r c and 0.04e/r c . In (y) and (z), curves 
from top to bottom correspond to f x — 0, 0.01e/r c , 0.02e/r c 
and 0.04e/r c . The horizontal dashed line corresponds to the 
value of Rg/yfi of a bulk solution at equilibrium. The vertical 
dot-dashed lines indicate the positions of the walls. 




0' i '- -L-^ -L 

-10 -5 5 10 

z[r c ] 




-10 -5 5 10 

z[r c ] 




FIG. 11: Distribution of the three components of chains ra- 
dius of gyration along the slit for the case of f x = 0.04e/r c and 
tp = 0.12. Curves from bottom to top correspond to TV = 5, 
20, 50, and 100, respectively. The vertical dot-dashed lines 
indicate the positions of the walls. 




FIG. 12: (a) Polymer density profile across the slit for the 
case of tp = 0.06 and TV = 100 for varying values of the driving 
force: f x =0 (squares), f x = 0.02e/Y c (diamonds), and f x = 
0.04e/r c (pluses), (b) Polymer density profile for the case of 
ip — 0.06, f x — 0.04e/V c , and chain length TV = 5 (circles), 
TV = 20 (diamonds), TV = 50 (pluses), and TV = 100 (squares). 



while the distribution is flat for small chains (TV = 5), a 
migration towards the center-plane is clearly observed as 
the chain length is increased (TV = 20), as demonstrated 
by an increase in the thickness of the depletion layers near 
the walls. However, at the very center-plane, the density 
profile exhibits a local minimum. As the chain length is 
further increased, the density profile in the center-plane 
is further decreased, while the two off-center peaks are 
increased, the depletion layer near the walls ultimately 
narrows, but slightly, as the chain length is increased. 
This rich behavior of the polymer density profile is qual- 
itatively in agreement with that reported by Usta et at 
using a lattice Boltzmann simulation 0], although in 
their case the size of the depletion layer in the weak con- 
finement regime seem to be strongly affected by the flow 
rate. 

The presence of a dip in the center-plane of the slit in 
Poiseuille flow implies a migration away from the very 
center of the slit. However, this can be accompanied by 
an increase in the depletion layer at the walls. Therefore 
an overall migration of the polymers cannot be merely 
inferred by merely looking at the depletion layer. We 
propose that the cross-stream migration issue can bet- 
ter be investigated through the second moment of the 
polymer volume fraction profile, 



H 



—H 



H 



dzz p p {z)l \ dzpp(z) 



-H 



(16) 



This is shown as a function of chain length, TV, for all 
systems investigated in Fig. 13. A small (z 2 ) implies a 
general migration from the walls. This figure shows that 
in the absence of flow (f x — 0), (z 2 ) decreases as TV in- 
creases. This is due to the increase in the depletion layer 
as TV is increased. In the presence of flow, (z 2 ) decreases 
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FIG. 13: The second moment as calculated from Eq. (16) vs. 
chain length. Dashed lines correspond to f x = 0, thin solid 
lines correspond to f x = 0.02e/r c ,and thick solid lines corre- 
spond to fx = 0.04e/r c . Black, blue and red lines correspond 
to ip — 0.06, 0.12, and 0.24, respectively. 

with increasing N for small values of if. However, for 
relatively large values of polymer volume fraction, (z 2 ) 
decreases with increasing N for small values of N, but 
for large values of N, (z 2 ) increases with increasing N. 
Furthermore, Fig. 13 shows that for small values of N, 
the level of migration increases with increasing the driv- 
ing force, f x . However, for large values of N, more mi- 
gration towards the wall is observed as f x is increased. 

V. SUMMARY AND CONCLUSIONS 



for non-Newtonian fluids with a power less than one. We 
also found that while the polymer chains are extended 
along the channel, the distribution of the three com- 
ponents of the radius of gyration is highly non-uniform 
across the channels. 

The density profile of the polymers in the presence of 
pressure driven flow is also found to be non-uniform, ex- 
hibiting a local minimum in the center-plane and two 
off-center symmetric peaks. We found either a migra- 
tion from or towards the wall depending on chain length, 
volume fraction and driving force. For relatively short 
chains and small polymer volume fractions, a migration 
away from the walls is observed. In this case, the level 
of migration is amplified as the chain length or driving 
force is increased. However, for longer chains a migration 
towards the wall is observed with the level of migration 
is increased as the chain length or the driving force is 
increased. For relatively concentrated polymer solutions, 
the level of migration is decreased. 

It was argued earlier by Groot and Warren j3^ that 
DPD suffers from a low Schmidt number, Sc = 77/ pDo ~ 
1, where r\ is the viscosity, p is the fluid density, and 
Do is the diffusivity of the solvent particles. In contrast 
a typical fluid has a Schmidt number of order 1000. A 
low Schmidt number implies that the particles diffuse 
as fast as momentum. Consequently, such an argument 
would lead us to conclude that the dynamics of polymer 
chains from DPD in a dilute solution will not follow the 
Zimm model. However, a recent extensive DPD study by 
us has shown that polymer chains in dilute solution cor- 
rectly obey the Zimm model. Therefore, within DPD, the 
hydrodynamic interaction is well developed within the 
time scale of polymer motion. As noted recently by Pe- 
ters [3^ , the Schmidt number in a coarse-grained model 
is in fact ill-defined, since in the Schmidt number, Dq 
corresponds to the diffusion coefficient of single particle, 
not coarse-grained fluid elements. Our current findings 
provide further testimony to the effectiveness of the dissi- 
pative particle dynamics approach in simulating polymer 
flow in a variety of situations including in nano-channels. 



We presented results obtained from a systematic com- 
putational study of flow of polymer solutions in nanoscale 
slit pores under the action of pressure gradient, using the 
dissipative particle dynamics approach. We found that 
while the flow of a simple fluid follows the Poiseuille- 
Hagen law, the flow of a polymer solution is characterized 
by a non-quadratic velocity profile and a lack of scaling 
with respect to both driving force and slit thickness that 
fits relatively well the generalized Poiseuille-Hagen law 
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